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> • ABSTRACT 

! The Centaurs are a population of small, planet-crossing objects in the outer 

solar system. They are dynamically short-lived and represent the transition pop- 
ulation between the Kuiper belt and the Jupiter family short-period comets. 
|"t"| ! Dynamical models and observations of the physical properties of the Centaurs 

! indicate that they may have multiple source populations in the trans-Neptunian 

region. It has been suggested that the inclination distribution of the Centaurs 
may be useful in distinguishing amongst these source regions. The Centaurs, 
I however, undergo many close encounters with the giant planets during their or- 

bital evolution; here we show that these encounters can substantially determine 
the inclination distribution of the Centaurs. Almost any plausible initial inclina- 
tion distribution of a Kuiper belt source results in Centaurs having inclinations 
peaked near 10 — 20°. Our studies also find that the Kuiper belt is an extremely 
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; unlikely source of the retrograde Centaur that has been observed. 

<N 

1. Introduction 

The group of minor planets known as Centaurs represent the dynamical link between 
the reservoirs of icy objects in outer solar system and the short-period Jupiter family comets 
in the inner solar system. Because the orbits of the Centaurs c r oss th ose of the outer planets, 



Centaurs are dynamically short-lived. iTiscareno fc Malhotral (120031 ) numerically integrated 



the orbits of 53 observed Centaurs and found that their median dynamical lifetime was only 
9 Myr. This means that there must exist a long-lived source for this transient population. 
Attempts to identify the Centaur sources have largely focused on dynamical models of the 
potential source regions. Early studies to identify a source region for the short period comets 
suggested that the long period comets (originating from the nearly isotropic distant Oort 
Cloud) could be captured into short period orbits by means of planetary perturbations, and 
that the higher capture probabiliti es at low inclina tions could expla i n the prog rade, low incli- 



nations of the short period comets (jEverhartlll972l ) . But iFernandea (119801 ) and lDuncan et al. 
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( 119881 ) showed that Everhart's mechanism is a very inefficient way to produce Centaurs and 
short period comets from Oort Cloud comets, suggesting that a more proximate, lower- 
inclination population just beyond Neptune would be a more likely source region. With the 
subsequent discovery of the Kuiper Belt and its complex dynamical structure, recent studies 
have focused on the dynamical subclasse s of the Kuiper belt as the sources that resupply the 
Cen taurs and the Jupiter fami ly comets: Levison Sz Duncan ( 1997 ). Volk fc Malhotra (I 2OO8 ) 
andlDi Sisto fc Bruninil ( 120071 ) modeled the scattered disk source; iDi Sisto et al.l (120101 ) and 
Morbidellil dl997h modeled the Plutinos (objects in the 3:2 mean motion resonance with Nep- 
tune) as a source of the short period comets, and the Trojan populations of Neptune and/or 
Ju piter (objects in the 1:1 m ean motion resonance with Neptu ne or Jupiter) were studied 
by iHorner fc Lykawkal (120101 ) and IHorner fc Wvn Evans! (12006). The Oort cloud has also 



been revisited as a possible source (fEmePyanenko et~ 2005 : Brasser et al. 2012 ). These 



models have identified dynamical pathways as well as the flux per unit source population 
from each source, but it remains unclear which, if any, source population dominates the flux 
of new Centaurs. We do not yet have strong enough observational estimates of the number 
of small, comet-sized (1-10 km) bodies in each of the proposed source regions to be able 
to esti mate the absolute number of new Centaurs and/or Jupiter family comets from each 



regio n ( iBernstein et al.l 12004c iDuncan Sz Levisonl 119971 : IVolk fc Malhotral l2008t iFraser et al. 
2010h . 



Another approach is to examine the dynamical and/or physical properties of the ob- 
served Centaurs and compare them to the same properties in the source region. Because 
the Centaurs have orbits that are much more favorable for observations than the orbits of 
trans- Neptunian obj ects (TNOs) , they have been the subject of several spectrosc opic and 
photometric studies (ITegler et al.l (120081 ) and references therein. iBarucci et al.l (120 1 lh ) . Their 
physical properties such as color and spectra l features are broa dly similar to those seen for 
the largest (and therefore brightest) TNOs (ITegler et al.l 120081 1 . supporting the dynamical 
connection between those populations. There are indications tha t the Centaurs hay e a bi- 
modal color distribution that might hint at multiple source regions ( ITegler et al.ll2008l ). More 
recently it has been suggested that the bimodal color distribution seen for the Centaurs is 
size dependent and can also be seen in small KBOs ( iPeixinho et al.ll2012l ). These obser- 
vations give us some clues about the origins of the Centaurs, but there are not yet large 
enough sample sizes to draw any firm conclusions. Furthermore, due to the inherent faint- 
ness of TNOs as observed from Earth, it is unclear that we will be able to observe sufficient 
numbers of smaller, Centaur-sized TNOs in the near future to improve this state. 

The dynamical properties of the Centaurs might be a possible way of identifying the 
source s. T he dynamics of the observed Centaurs have been modeled by lTiscareno fc Malhotra 
(120031 ) and iBailey fc Malhotral (120091 ): the dynamics of the known Centaurs is found to be 
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dominated by chaotic diffusion and the frequent close encounters with the outer planets, 
which cause the sem imajor axes and ec centricities of the Centaurs to evolve rapidly. It has 
been suggested (e.g. iGulbis et all (120 lOh ) that the orbital inclinations of the Centaurs might 
preserve some memory of their source regions. There are currently sufficient observations 
of TNOs to identify several dynamical subclasses and to calculate inclination distributions 
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observed inclination distribution of these objects is sho wn in Fig. [T] along with the debiased 



inclination distribution from the Deep Ecliptic Survey ( IGulbis et al.ll2010l ). If the Centaurs' 
inclination distribution does retain some memory of the source region's inclination distribu- 
tion, that could be an important means for linking the Centaurs to a particular dynamical 
subclass of the Kuiper belt or to the Oort cloud. 

In this work, we investigate the extent to which inclinations of Centaurs may preserve 
the inclination distribution of their sources. We present the results of numerical simulations 
of hypothetical Centaurs generated from the Kuiper belt subclasses. We follow the orbital 
evolution of these hypothetical Centaurs until they either transition onto inner solar system 
orbits or they are ejected from the giant planet region. During this evolution, our simulated 
Centaurs experience many close encounters with the four outer planets, each of which induces 
some change in the orbital elements of the Centaurs. We use these simulations to investigate 
the average outcomes of planetary encounters, including the average change in a Centaur's 
inclination as a function of close encounter distance. Using these average encounter outcomes, 
along with statistics about close encounter frequency as a function of encounter distance 
and the dynamical lifetimes of Centaurs, we compute the expected dynamically evolved 
inclination distribution of Centaurs as a function of the initial inclination distribution of 
the source population. We then discuss the implications of these results for identifying the 
dominant Centaur sourcefs). 



2. Numerical Simulations 

To explore the evolution of inclinations in the Centaur region, we performed a numerical 
integration of test particles entering the planet-crossing region from a trans-Neptunian source 
region. We chose source region parameters based on the classical Kuiper belt (CKB) and the 
scattered disk (SD) because these dynamical classes are widely regarded as the most likely 
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tive of either the CKB or the SD, but is rather an exploration of the parameter space of 
these populations. To that end, our test particles uniformly cover the following ranges in 
semimajor axis (a), perihelion distance (q), mean anomaly (M), argument of perihelion (oj), 
and longitude of ascending node (fi): 



• 40 AU < a < 80 AU 

• 30 AU < q < 32 AU 

• < M, fi, w < 2vr. 



To explore a large range in inclination, the test particles were divided into 24 inclination 
bins (200 test particles per bin) with the following bin widths: 



• 1° for 0° < i < 10° 

• 2° for 10° < i < 30 

• 5° for 30° < i < 50 



and inclinations were assigned uniformly within each bin. The choice to limit perihelion 
distances to the 30 — 32 AU range was made to manage computational time. Test particles 
with these q values will either enter the Centaur population or be ejected from the solar 
system early in the simulation, minimizing the necessary total length of the simulation. We 
also limit the initial semimajor axes to < 80 AU to conserve computational time because test 
particles starting at larger values of a will take longer to evolve onto Centaur-like orbits. We 
show later in this section that the dynamical properties we are interested in do not depend 
on this limit in a. 

The integration was performed using the swift_rmvs3 code in the SWIFT software pack- 
agdE which is capable of integrating test particles through arbitrarily close encounters with 
the massive planets. For each simulation we include the four outer planets and the sun as 
massive bodies (the sun's mass is augmented with the mass of the terrestrial planets), and the 
hypothetical Centaurs are included as massless test particles. We use a timestep of one year, 
and we remove test particles if they achieve a heliocentric distance inside Jupiter's orbit or 
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beyond 1000 AU, or if they impact one of the giant planets. The inner boundary of the simu- 
lation was chosen to conserve computational time; following test particles into the inner solar 
system would require the inclusion of the terrestrial planets and the use of a much shorter 
timestep. This choice means that we are not including in our fictitious Centaur population 
those objects that enter the inner solar system, become Jupiter family comets, then evolve 
back onto Centaur orbits. Previous work has shown that once test particles evolve to Jupiter 
crossing orbits they are either ejected from the solar system or they become short-period 
comets which are then ejected from the solar system on very short timescales (~ 10 5 yrs) 
f lTiscareno fc Malhotr albood : bi Sisto fc Bruninil 120071 ; bailev fc Malhotrall2009f ) ; this means 
that for time-weighted distributions in the Centaur region, our inner boundary will have al- 
most no effect on the results of our simulations. The simulation was run for 300 Myr because 
this is many times longer than the expected dynamical lifetime (~ 10 Myr) for most Cen- 
taurs a nd longer than even t he ve r y long lived Centaurs simul ated by iTiscareno fc Malhotra 
f l2003h : bi Sisto & Bruninil J2007h : iBailev k Malhotra! feood ). In the analysis of our simu- 
lations, we are particularly interested in the dynamical lifetimes of the Centaurs produced 
from our source population, as well as the close encounters they have with the four outer 
planets and the evolution of the inclinations of the test particles. In Section [3] we discuss 
how we use the results of this simulation to map arbitrary source inclination distributions to 
their Centaur inclination distribution without the need for additional numerical simulations. 



2.1. Planetary encounter statistics 

The output for our simulation includes the orbital histories of all the test particles as 
well as the detailed outcomes of more than 1.8 million planetary encounters with Neptune, 
Uranus, Saturn and Jupiter. A planetary encounter is defined as approaching within one 
Hill radius of a planet. During its lifetime as a Centaur, each test particle has an average 
of 420 planetary encounters before it either enters the inner solar system or is ejected from 
the solar system to distances greater than 1000 AU. Fig. [2] shows the distribution of closest 
approach distances, r ce . The distances are scaled to the encountered planet's Hill radius, 
given by: 

R H = a p ( , (1) 

where a p is the planet's semimajor axis, M p is its mass, and M Q is the mass of the sun. 
As Fig. [2] shows, 50% of the encounters in the simulation occur at separations larger than 
0.7Rh- The distribution of closest approach distances follows a power law ~ (i?///r ce ) 2 . 
This is consistent with expectations from geometrical arguments that neglect gravitational 
focusing due to the planet. 
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We also examine the change in a test particle's heliocentric velocity, Av = \v) — v~l\, 
imparted by each planetary encounter. In the limit of test particles with nearly circular, 
coplanar orbits very close to a planet, an analytic estimate of Av is as follows. Let the 
particle's semimajor axis be a = a p (l+x); the closest approach distance is given by r ce = xa p . 
Let us define the encounter time, t enc , to be the time it takes to move an azimuthal distance 
2x. Then 

2x 



v enc 



n — n p 



(2) 



where n and n p are the mean motions of the test particle and planet respectively. Assuming 
that x << 1, we can estimate the test particle's mean motion as 

n = n p (l + x)~ 3/2 w --n p x (3) 

which gives us t enc = The velocity kick from the planetary encounter can then be 

estimated as 

Av = 



Grrip 



{xa p ) 2 



W = ^(Gm ) 1/2 afr c - e 2 (4) 
3 m ' 

Scaling the encounter distance to the planet's Hill radius, we expect nearly circular test 
particles to experience a Av that scales with planet mass and semimajor axis in the following 

way 



The orbits of the Centaurs are neither circular nor coplanar, but this serves as a useful 
comparison for the simulation results. The left panel of Fig. [3] shows the average absolute 
value of Av as a function of r ce /R H ; we have divided the encounters by planet, and the 
error bars on the plot show the rms variations in Av for each bin in r ce /Rn- The majority 
of encounters have only very small effects on the test particles' orbits because they occur 
at relatively large separations. The average Av as a function of encounter distance can be 
fit to a power law in the scaled closest approach distance, the mass of the planet, and the 
semimajor axis of the planet: 

The fit parameters were determined using the simulated close encounters with Saturn, 
Uranus, and Neptune; Jupiter encounters were excluded from the fits because the inner 
boundary of the simulation does not allow us to follow a sufficient number of test particles 
through Jupiter encounters. However, Fig. [3] shows the averaged results of the encounters 
for all four planets along with the above power law; the results for Jupiter are consistent 
with Eq. El 
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The simulation results are consistent with the analytic estimate (Eq. [5]) for the scaling 
of Av with r ce , but not for the scaling in planet mass and semimajor axis. The mean \Av\ 
appears to be nearly insensitive to the planet mass. These differences are partially due to 
the substantial eccentricities and inclinations of the Centaurs in our simulation; this affects 
their encounter velocities and encounter times. The range of planet masses is quite small: 
we have only three planet masses in our fit, Saturn, Uranus and Neptune, and the latter 
two have very similar masses, which means that the planet mass dependence of our best-fit, 
Eq. is not very well determined. The planet mass and semimajor axis scalings could also 
be dependent on the architecture of the giant planets in our solar system; a different planet 
configuration could yield a different scaling. 

The change in a test particle's inclination due to planetary encounters is determined 
by the fraction of Aw that is in the z direction. We can again use the simulation results to 
determine the average change in inclination, Ai as a function of the planet mass, semimajor 
axis, and the closest approach distance: 

The change in the exponent for the a p /AU term as compared to Eq. [6] is due to the fact that 
Ai oc Av z /v, and the orbital velocity of the test particles during a planetary encounter is 
roughly proportional a~ 0,5 . 



2.2. Dynamical lifetimes and inclination changes 



We also exa mine the dynamical lifetimes of the Centaurs in our simulation. As noted in 
previ ous studies (IDi Sisto fc Bruninil l2007t iTiscareno fc Malhotral l2003t iBailey fc Malhotra 
20091 ). the Centaurs show a v e ry wi de range of lifetimes, from only ~ 1 Myr to more than 
100 Myr. IDi Sisto fc Bruninil (120071 ) noted that a Centaur's lifetime depends on its initial 
inclination (with larger inclinations correlated with longer lifetimes). Fig. H] shows the median 
dynamical lifetimes, r, of our simulated Centaurs as a function of their initial inclination in 
the source region; the binning scheme is the same as that for the initial conditions. Only 
test particles that entered the Centaur population (5 < q < 30 AU), evolved, and were 
subsequently removed from the simulation (at the inner or outer boundaries, or after a 
collision with a planet) are included in the calculation of the dynamical lifetimes. Even with 
our initial, near-Neptune perihelion distances, some test particles don't enter the Centaur 
population until late in the simulation; because these test particles are still evolving as 
Centaurs when the simulation was stopped, we do not have accurate values of r for them, so 
they are excluded. The time spent as a Centaur increases with initial inclination, with the 
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highest inclination test particles spending about 4 to 5 times longer as Centaurs than the 



l owest inclination test particles. This is consistent with previous work by lDi Sisto fc Brunini 



( 120071 ) (see their Fig. 5), although the values of r are slightly different because they use the 
mean value of r while we use the median. We prefer to use the median value for the 
dynamical lifetimes because the distribution of r for any bin in initial inclination has a very 
long tail, which tends to skew the mean toward larger r. Fig. H] also shows the distribution 
of r for three different ranges in initial inclination. The increase in r with increasing initial 
i is simply due to a decrease in the rate of planetary encounters because high-2 test particles 
spend little time near the plane of the planets. The number of planetary encounters required 
to either transfer a test particle onto an orbit that enters the inner solar system or eject it 
does not depend on initial inclination (for the range of inclinations we tested), so the lower 
frequency of encounters corresponds to longer dynamical lifetimes. 

The net result of all these planetary encounters can be seen in Fig. [5] which shows the 
average inclinations of the Centaurs in our simulation as a function of their initial inclina- 
tions. Test particles that start with inclinations below ~ 10° have their inclinations nearly 
uniformly raised to ~ 10°. At higher initial inclinations, the planetary encounters have a 
net effect of slightly lowering the average inclinations. To see how dependent these results 
are on the semimajor axis distribution of our source region, we performed an additional 
simulation with 800 test particles spread uniformly in semimajor axis from 100 — 500 AU 
with the same inclination and perihelion distance distribution as outlined above. We cal- 
culated the average inclinations and median dynamical lifetimes for the resulting Centaurs 
and performed a Spearman rank correlation test to see how these two parameters depend 
on the initial semimajor axis. This test produces a number in the range [—1, 1], with results 
near meaning that the two variables being compared are not correlated (changes in one 
variable do not predictably result in changes to the other). Comparing initial semimajor 
axis with Centaur lifetimes and with the average change in inclination during that lifetime, 
we find Spearman rank correlation coefficients of ~ 10~ 2 , consistent with no correlation be- 
tween initial a and either dynamical lifetime or Ai. This means that the averaged results 
from our initial 2400 test particle simulation should be applicable to trans-Neptunian source 
regions with semimajor axis distributions that differ from our uniform 40 < a < 80 AU test 
simulation. 



3. Mapping the Source Inclination Distribution to the Centaur Region 



The results of the simulation in Section [2] tell us, on average, how the inclination and 
dynamical lifetime of a Centaur is related to its initial inclination. As discussed in Section (2j 
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these two key parameters are not much affected by the initial semimajor axis distribution 
of the test particles, at least within the range of semimajor axes typical for the Kuiper belt 
and the scattered disk. This means that without worrying about the distributions of the 
other orbital elements, we can take an arbitrary source region inclination distribution and 
calculate what the ^-distribution would be for the resulting Centaurs. We divide the range 
of initial inclinations into bins (using the same scheme as in Section [2]), and we model the 
inclinations of the Centaurs in each bin as a Gaussian multiplied by sin(i): 



g(i, fi,a) = C sin(i) exp 



2a 2 



where \i and a are fit parameters determined from the simulation results and C is a nor- 
malization constant; the values of these parameters are listed in Table [TJ Using these fits to 
the simulation results and the median dynamical lifetimes for each inclination b m , T n , as a 
weighting factor (also listed in Tabled]), we calculate the time-averaged Centaur inclination 
distribution as follows: 

ribv 



"-bins r ^ rl-n+1 

f centaur ji) = C ^ ] g n {% , /!, 0")r n — — / f s 

„_i L l^n+l ~ In) J i_ 



[i)di 



(9) 



n=l 

where i n is the lower bound on the initial inclinations in bin n, and C is a normalization factor. 
Fig. [6] shows the results of this mapping procedure for two trial source region inclination 
distributions; both source distributions take the form of a Gaussian multiplied by sin(z) 
with the left panel representing a very low-? source region and the right panel representing 
a higher-z source region. For a very low-2 source region, the resulting Centaur ^-distribution 
has a much larger spread and a larger mean i than the source; the planetary encounters 
systematically raise the inclinations of the Centaurs. Centaurs from a high-z source region 
will have an inclination distribution similar to the source region, but skewed slightly toward 
lower i. 



Comparing the Mapping Procedure to the Results of Numerical 

Simulations 



To test the accuracy of the mapping procedure, we performed two additional numerical 
simulations: one of a CKB-like source region and one of a SD-like source region. The CKB- 
like source consisted of 2400 test particles with semimajor axe s spread evenl y over the range 
40 AU < a < 50 AU with an inclination distribution given by iBrownl (120011 ) 's debiased best 
fit to the observed classical belt objects: 



smm 



0.93 exp 



2a( 



0.07 exp 



1a\j 



at = 1.4°, a 2 = 17.0 C 



(10) 
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As with the initial simulation in Section |2l the perihelion distances were evenly spread in 
the range 30AU < q < 32AU. This distribution of test particles is by no means meant to 
be a representation of the true CKB orbital distribution; instead it is merely a test of the 
applicability if our initial simulation results to a CKB-like orbital parameters and a realistic 
CKB inclination distribution. 

The SD-like source region consists of 2400 test particles spread over a semimajor axis 
range 50 < a < 250, perihelion range 30AU < q < 32AU, and an inclination distribution 
given by: 

/(z) = sin(z)exp(^,a = 17°, (11) 
which is the debiased best fit inclination distribution for the observed set of scattered disk 



objects (IVolk fc Malhotrall2008l ). Again, this SD-like source population is merely a test of 
the mapping procedure, not an accurate portrayal of the scattered disk distribution. Both 
sets of source region test particles were integrated using the same simulation parameters 
described in Section [2j Fig. [7] shows the initial ^-distributions for these simulations, the Cen- 
taur distribution predicted by the mapping procedure, and the time-averaged i-distribution 
for the Centaurs from the simulation. The mapping procedure is quite successful in both 
cases, giving confidence that the mapping procedure can be used in lieu of detailed simula- 
tions to predict the Centaur ^-distribution that would result from an arbitrary source region 
distribution. 



5. Results 

5.1. Inclination as a discriminator of Centaur sources 

It is apparent from the simulations that repeated encounters with the planets will signif- 
icantly alter the inclination distribution of the Centaurs. If there are any peaks in the source 
region's inclination distribution at i < 10°, these will generally not be preserved as the incli- 
nations of the new Centaurs are raised by planetary encounters. This is particularly relevant 
for Centaurs that might originate in the classical Kuiper belt, which has been shown to have a 
double peaked inclination distribution. Even if Centaurs leak out of the CKB independently 
of their original inclinations, the double peaked nature of the inclination distribution will 
not persist in the Centaur region. The reality is even less favor able for preserving the double 



peaked distribution, however. In our model of the CKB from IVolk fc Malhotral (1201 if ), we 
found that CKBOs from the high-i portion of the inclination distribution were ~ 4 times 
more likely to evolve onto Neptune crossing orbits than CKBOs from the \ow-i peak. The 
result is that trying to use the Centaur i-distribution to distinguish whether the majority 
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of Centaurs originate in the scattered disk or in the classical belt will be nearly impossible 
because the high-z CKBOs and the SDOs have very similar inclination distributions; in the 
case of our test simulations in Section HJ the Centaurs from the SD-like source region are 
nearly identical to those from the CKBO-like source region. 

The best chance to be able to identify source regions using the Centaur inclination 
distribution is if the proposed source regions have well separated high-inclina tion features. 



as the se are better preserved than low-* features. Taking the i-distributions from lGulbis et al. 



(120101 ) for several subclasses of TNOs as examples, the difference between the ^-distribution 
of the scattered objects (fi = 19°, a = 7°) and the resonant objects (/i = 6°, a = 9°) might 
be large enough to be observationally distinguishable in the Centaur region. Fig. [8] shows 
these two different source region inclination distributions mapped to their resulting Centaur 
inclination distributions using the method described in Section [3J Referring to Fig. [TJ it 
is clear that there are far too few debiased Centaur observations to be able to make this 
distinction at the present time. We estimate that the two distributions in Fig. [8] would be 
observationally distinguishable with approximately 100 unbiased Centaur detections with 
< i < 30°; this estimate is based on binning the distributions into 5° bins, assuming 
Poisson noise statistics, and requiring that the la error bars in the resulting histograms 
overlap in fewer than half the inclination bins. 



5.2. Retrograde and very high inclination Centaurs 



We mined our simulation results to identify cases where individual test particles evolved 
to very high inclinations in the Centaur region. This is motivated by the di scoveries of 



sever al Centaurs with very large inclin ations, such as 2002 X U93 with i ~ 78° ( lElliot et al. 
2005h and 2008 KV42, with i ~ 104° rtdadman etaDl2009h . We examined our simulation 



results to see if test particles with lower i Q (typical of TNOs) could be excited to such large 
inclinations. Our simulated Centaurs almost all started with i < 50°, and although some 
achieved inclinations as large as 80 — 90° within the Centaur region, they did so only very 
briefly (timescales < 1 Myr) before being ejected on hyperbolic orbits. These also tended to 
have very large semimajor axes, a > 1000AU, when they were at large inclinations, which 
is inconsistent with the value a ~ 40 AU found for 2008 KV42. We numerically integrated 
nearly 10,000 TNO-like test particles in the simulations presented in this work, and after 
following them all through many encounters with the planets, we find no Centaurs with 
orbits similar to t he two previously m entioned high-z Centaurs. This result is consistent 
with the results of iBrasser et al.l (120 121 ) who find the probability of evolving from a Kuiper 
belt-like orbit to a nearly retrograde Centaur orbit is ~ 10 -5 . Therefore it appears unlikely 
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that these relatively stable (stable on ~ 10 Myr timescales) observed high-i Centaurs can be 
satisfactorily explained merely by typical orbital evolution of TNOs due to close encounters 
with the planets. This can be understood quantitatively given that planetary encounters 
roughly preserve the Tisserand parameter, 



T 



+ 2 



1/2 



cosz. 



(12) 



Large changes in cos i (or a change in the sign of cos i as in the case of prograde to retrograde 
orbits) will tend to be associated with changes in a, possibly removing it from the C entaur 
region. An alternate, very high-z source region, as discussed by lGladman et al.l ( 120091 ). might 
be more likely. 

We performed an additional numerical simulation of test particles with very large incli- 
nations, 50 < i < 90, and semimajor axes in the range 50 < a < 250 AU. In the simulation 
we found that Centaurs with io >~ 70° were able to evolve to i ~ 100°. As in our other 
simulations, most of these high values of i were achieved only briefly while at very large 
values of a. But we did find one case where a test particle with i = 75° spent ~ 10 Myr 
with an orbit very similar to that of 2008 KV42 (a ~ 25 AU, i ~ 105°, q ~ 20 AU); this case 
is shown in Fig. [91 This simulation demonstrates that there does exist a dynamical path- 
way from very high-z, scattered disk-like orbits to quasi-stable retrograde Centaur orbits. A 
robust estimate of the likelihood of this evolution is beyond the scope of the present work, 
but based on the 237 test particles that completed their evolution as Centaurs during the 
simulation, the probability of evolving from 50 < i < 90° onto an orbit similar to 2008 KB42 
is ~ 4 x 10~ 3 . 

A potentially more likely source regio n for these pa r ticula rly high-inclination Centaurs 
is the Oort cloud, as recently suggested by lBrasser et al.l ( 120121 ). We did a test simulation of 
Centaurs coming from an isotropic inclination distribution with semimajor axes in the range 
50 — 250 AU (consistent with the semimajor axis range of test particle s that have entered t he 
planet crossing region from the Oort cloud, see for example Fig. 2 in iBrasser et al.l ( 120121 )). 
We find that the Centaurs originating on retrograde orbits have median lifetimes of ~ 250 
Myr and have average inclinations very similar to their source region inclination. This sug- 
gests that any contribution to the Centaur population from a source region with retrograde 
inclinations should be long-lived compared to Centaurs originating from a prograde source 
region. 
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6. Summary 

In summary, we have performed numerical simulations of the dynamical evolution of 
Centaurs originating from the Kuiper belt and the scattered disk and have found that: 

• Multiple planetary encounters can significantly alter the inclination distribution of 
Centaurs relative to their source region's inclination distribution. 

— Test particles with < 10° tend to have their inclinations raised to values near 
10°. 

— Test particles with z'o > 10° tend to have their inclinations lowered by a small 
amount that still preserves the correlation between average inclination and initial 
inclination. 

— High-inclination features of the source region are better preserved than low- 
inclination features. 

• The median dynamical lifetimes of Centaurs increase nearly linearly with initial incli- 
nation. 

• The average absolute change in velocity of Centaurs due to planetary encounters can 
be represented as a power law, |Au| = 6.24 x 10 3 (ms -1 )(^-)°- 18 (|| F )- 2 - 1 (^;)- 2 . 

• Planetary encounters are not an efficient way to produce extremely high-inclination 
Centaurs. However, there is a dynamical pathway from the scattered disk at z'o ~ 70° 
to retrograde Centaur orbits that are stable for several million years. 

This research was supported by grant no. NNX08AQ65G from NASA's Outer Planets Re- 
search program. 
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Fig. 1. — The inclination distribution of the observed Centaurs (taken from the Minor Planet 
Center website) and the debiased inclination distribution for the Centaurs from the Deep 
Ecliptic Survey (IGulbis et al.ll2010l ). 
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Fig. 2. — Cumulative fraction of encounters (solid line) occurring at closest approach dis- 
tances less than r ce jRn. The dashed lines indicate the approach distance outside of which 
50% of encounters occur. 
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Fig. 3. — Average \Av\ per encounter (left panel) and average \Ai\ per encounter (right 
panel) vs the scaled closest approach distance for the simulation results (solid lines) and the 
power law fits (dashed lines, Eq. [6] and Eq. [7]). The error bars in the left panel are the rms 
variations in \Av\. 
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Fig. 4. — Left panel: median dynamical lifetimes as a Centaur vs. initial inclination in the 
source region. The bin sizes vary such that there are roughly equal number os test particles 
contributing to each bin. Right panel: distribution of Centaur lifetimes for three different 
initial inclination ranges. 
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Fig. 5. — Average inclination as a Centaur vs initial inclination in the source region. The 
error bars are the rms variation in i within each initial inclination bin, and the dashed line 
show a one to one correlation for reference. 
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Fig. 6. — Inclination distributions for a trial source population with an inclination distri- 
bution given by Eq. [8] with /i = 2° and a = 2° (left panel) and fi = 19° and o = 7° (right 
panel) compared to the Centaur distributions given by the mapping procedure (Eq. [9]). 
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Fig. 7. — Initial inclination distributions, the predicted Centaur distribution calculated from 
Eq. |Hl and the time averaged Centaur distribution from a simulation for a SD-like source 
population with an inclination distribution given by Eq. (TT] (left panel) and a CKB-like 
source population with an inclination distribution given by Eq. [TU] (right panel). 




Fig. 8. — Predicted Centaur distributions calcu lated from Eq. |9] for initial source region 
inclination distributions from iGulbis et al.l (120101 ) for the scattered disk (solid line) and the 
resonant Kuiper belt (dashed line). 
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Fig. 9. — Evolution of a test particle with iq = 75° onto a quasi-stable orbit (shaded region) 
similar to the observed retrograde Centaur 2008 KV42. The values of a, e, and % for 2008 
KV42 are shown by the dashed lines while the evolution of the test particle is shown in solid 
lines. 



Table 1. Fit parameters for Eq. [8] applied to the simulation results for each bin in initial 
inclination, iq (see Section [2]), along with the median dynamical lifetime, r, of Centaurs 

from each bin in iq. 
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